function y = func_J_m2(q, m, CONSTS)

    check_asympt = false;
    k0 = CONSTS.k0;
    a = CONSTS.a;
    p = func_p_from_q(q, CONSTS);
    q2 = func_q2_from_p(p, CONSTS);
    beta_2 = func_beta_2(p, CONSTS);

    y = (besselj(m+1, k0*a*q2)./((k0*a*q2).*besselj(m, k0*a*q2))) - ...
           ((beta_2.*m)./((k0*a*q2).^2));
       
    if(check_asympt)  % asymptotics is held for q > 1000 rigorously 
        eps = CONSTS.eps;
        eta = CONSTS.eta;
        g = CONSTS.g;
        beta_2_s = 1+(g/(eta-eps));
        figure; plot(q, real(beta_2), 'b.-', q, real(beta_2_s), 'r-'); title('Re(\beta_2)');
        figure; plot(q, imag(beta_2), 'b.-', q, imag(beta_2_s), 'r-'); title('Im(\beta_2)');
    end

end
